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We discuss how dynamical fermion computations may be made yet cheaper by using symplectic 
integrators that conserve energy much more accurately without decreasing the integration step 
size. We first explain why symplectic integrators exactly conserve a "shadow" Hamiltonian close 
to the desired one, and how this Hamiltonian may be computed in terms of Poisson brackets. We 
then discuss how classical mechanics may be implemented on Lie groups and derive the form of 
the Poisson brackets and force terms for some interesting integrators such as those making use of 
second derivatives of the action (Hessian or force gradient integrators). We hope that these will 
be seen to greatly improve energy conservation for only a small additional cost and that their use 
will significantly reduce the cost of dynamical fermion computations. 
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1. Symplectic Integrators 



We are interested in finding the classical trajectory in phase space of a system described by the 
Hamiltonian H{q,p) = T{p) +S{q) = ^p^ +S{q). The idea of a symplectic integrator is to write 
the time evolution operator as 



exp T 



dt 



■ exp T 



dp d dq d 
dt dp dt dq 



where the vector field 



H 



dH d dH d d d . 

dq dp dp dq dp dq 



Since the kinetic energy T is a function only of p and the potential energy 5 is a function only of 
q it follows that the action of e'^^ : f{q,p) ^ f{q,p — TS'{q)) and e'^^ : f{q,p) ^ f {q + zT' (p) , p) 
are just translations of the appropriate variable. 

We now make use of the Baker-Campbell-Hausdorff (BCH) formula, which tells us that the 
product of exponentials in any associative algebra can be written as ln(e^/^e^e^/^) — {A + B) = 

^ |[A, [A,B]] — 2[B, [A,B]]| H where all the terms on the right hand side are constructed out 

of commutators of A and B with known coefficients. We find that for a simple PQP symmetric 
integrator with step size 5z the evolution operator for a trajectory of length T may be written as 



t/St 



exp 
exp 



f + S)5T-^J[S,[S,f]]+2[f,[S,f]]]5z' + 



z(f + S-^A[S,[S,f]]+2[f,[S,f]]]5T' + &{5x') 



2. Shadow Hamiltonians 



For every symplectic integrator there is a shadow Hamiltonian H that is exactly conserved; 
this may be obtained by replacing the commutators [S, f] in the BCH expansion with the Poisson 
bracket {S,T} = — f^f^ For example our PQP integrator above exactly conserves the 

shadow Hamiltonian H = T + S-^, {{S, {5, T}} + 2{r, {5, T}}) St^ + . . .. 

We now make the simple observation that any symplectic integrator is constructed from the 
same Poisson brackets, and that these Poisson brackets are extensive quantities. We therefore 
propose to measure the average values of the Poisson brackets and then optimize the integrator 
(by adjusting the step sizes, order of the integration scheme, integrator parameters, number of 
pseudofermion fields, etc. ||^, ^) offline so as to minimize the cost. This is possible because the 
acceptance rate and instabilities are completely determined by 8H = H — H. 

As a very simple example consider the minimum norm PQPQP integrator 
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Flat Manifold 


General 


Symplectic 2-form 
Hamiltonian vector field 
Equations of motion 
Poisson bracket 


dp Adq 

^ _ dH d dH d 
dp dq dq dp 

L ' J dp dq dq dp 


O):dco = 
dH = ifjCO 

d 1 M 

Tt\o-^ 
= -(o{A,B) 



Table 1: Comparison of quantities in flat space and on a general manifold 



whose shadow Hamiltonian is 



6a^-6a + l 1 -6a \ ^ 2 



H = H+ {S. T}} + —^{T. {S, T}}j 5x' + ^(5t^). 

With only one degree of freedom a we cannot completely eliminate the coefficient of the 0{5x^) 

I! J J\\\ 

contribution, however, we may optimize this integrator by setting the parameter a = { + kj{i\jfY^- 
There have been alternative optimization strategies proposed: minimizing the L2 norm of coeffi- 
cients assuming | {S, {S, T}}\ = \ {T, {5, T}} | [Q], and setting the coefficient of one of the two Pois- 
son brackets to zero by choosing a = ^(l — ■^)org. However, these strategies clearly break down 
when optimizing higher order minimum norm integrators, i.e., for C?(5t^) integrators there are 6 
Poisson bracket contributions that must be considered (see Table ph. 



3. Hessian Integrators 

We now make another simple observation: consider again the PQPQP integrator, where we 
set a = g so that the {T, {5, T}} contribution is eliminated. The remaining leading order Pois- 
son bracket {5, {5, T}} depends only on q, which means that we can evaluate the integrator step 
^{s,{S-T}}8t^ explicitly (it is again just a shift of p). The force for this integrator step involves sec- 
ond derivatives of the action, and therefore they are called Hessian or force gradient integrators 
[Q, ^. By putting such an integration step into a multistep integrator we can eliminate all the lead- 
ing ^(5t^) terms in 8H. The advantage of such an integrator over that of Campostini ||] is that 
the coefficients of the next order terms are approximately two orders of magnitude smaller (see Ta- 
ble |3|). We want to stress that although eliminating the leading term must be best asymptotically as 
5t — > it might well not be the optimal solution in practice; the optimal solution may be obtained 
by minimizing 8H as discussed in §§. 

4. Beyond Scalar Field Theory 

We now have to construct the Poisson brackets and Hessian integrators for gauge fields, where 
the field variables are constrained to live on a group manifold. To do this we need to use some 
differential geometry. Table |l| summarizes the difference between the formulation on fiat space 
that we have discussed up to this point and that on general manifolds. 
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In order to construct a Hamiltonian system on a manifold we need not only a Hamiltonian 
function but also a fundamental closed 2-form CO. On a Lie group manifold this is most easily 
found using the globally defined Maurer-Cartan forms {6'} that are dual to the generators and 

are the structure constants of the group. We 



where c 



jk 



satisfy the relation dd' 

choose to define co = -JL- d'p' = l,i{d' Adp' - p'dd') = l,i{d' Adp' + ^p'c'j^di A d"). Using 
this fundamental 2-form we can define a Hamiltonian vector field A corresponding to any 0-form 



A through the relation dA 
this gives 



i^co, and in the natural coordinates 



dA 
dp 



■ei + 



jk 



dA 
Ipi 



■ei{A) 




on the contangent bundle 



(4.1) 



The classical trajectories Ot = {Qt,Pt) are then the integral curves of this vector field, dt =A{Ot). 



5. Putting It All Together 

Recalling that H = S + T v/e can compute the Hamiltonian vector fields corresponding to S 
and T using equation (4. 1 ), and from these we can evaluate the lowest-order Poisson bracket 

dS 



{S,T} = -(0{S,f) = -{d' Adp' + yc)kdj Ad''){S,f) = -p'eiiS) 

and the Hamiltonian vector field corresponding to it, 

'd{S,T} 



dp' 



-ei + 



jk 



dpi 



ReTr 



d_ 
Tp' 



dU 



PU 



-ei{S)ei 



-%p'ej{S)+phiej{S) 



dp' 



From this we can derive expressions for the third- and fifth-order Poisson brackets that are needed 
for symmetric symplectic integrators, and these are listed in Table Similarly, we can then evaluate 
the corresponding Hamiltonian vector fields for any Poisson brackets we wish to include in the 
integration (e.g., {5',{5,r}} for force gradient integrators). 

The explicit form of the shadow Hamiltonian for a variety of integrators is show in Table ^. 



6. Conclusions 

Our work in this area is still very preliminary, so far we have concentrated on developing these 
ideas. Future work shall focus on implementing and testing the performance of these integrators 
for dynamical fermion calculations. We expect that modest gains in performance can be expected 
through directly measuring the leading order Poisson brackets to optimize the minimum norm 
family of integrators. However, we hope that very significant performance improvements can be 
obtained from force gradient integrators. 
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{T,{S,T}} 
{S,{S,{S,{S,T}}}} 
{{S,T},{S,{S,T}}} 

{{S,T},{T,{S,T}}} 

{T,{S,{S,{S,T}}}} 
{T,{T,{S,{S,T}}}} 
{T,{T,{T,{S,T}}}} 



eiiS)ei{S) 
-p'p%ej{S) 


-2ei{S)ej{S)eiej{S) 
3c'ji^p'p^ej{S) [ekefXS) + efe^iS)] 

fnJ 



+P'P 



ek{S)ekeiej{S) - [ekei{S) + eiek{S)]eke j{S) 


2p'p^ [eie jek{S)ek{S) + eiek{S)e jek{S)] 
-p'p-'p''p^eiejekee{S) 



Table 2: Poisson brackets required for symmetric symplectic integrators. 
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Integrator 


Update steps 


Shadow Hamiltonian 


PQP 


^StS UtT ^StS 


T + S-^ ({5, {5, T}} + 2{r, {5, T}} 


) 


QPQ 




T + S+^[2{S,{S,T}} + {T,{S,T}} 


) 


PQPQP 

Sii 


X £3^^^ 


T + S+^{S,{S,T}} 


PQPQP 




T + S+^5t^{T,{S,T}} 


Campostrini 


exp(^+^f+^5Tf) 
X exp('^+^g^+^5T5) 


) 
) 
) 


T + S 

( -(40^+40^+48) {5, {5, {5, {5, r}}}} ^ 

+ (180^+240^+312) {{S, T], {S, {5, r}}} 

+(60^+80^+104) {{5,r},{r,{5,r}}} 

~r 34560 

+(-2oe^+8){r,{5,{5,{5,r}}}} 
+(2o^+32){r,{r,{5,{5,r}}}} 
V +(5^+8){r,{r,{r,{5,r}}}} / 


X exp(-^-/;^+^5Tf 
xexp( ^+^3^+'5t5 
X exp(-'^-/;^+^5Tf 
X exp('^+^g^+^5T5) 
X exp('^'^ff^'5Tt) 


Force 
Gradient 
#1 

M 


eiSrt ^ISrS ^^Srt 
X e 192 


_1_ St4 
~r 6635520 


T + S 

( 2259 {5, {5, {5, {5, r}}}} ^ 
+3024 {{5, r}, {5, {5, r}}} 
+768 {{5, r},{r, {5, r}}} 

+5616{r,{5,{5,{5,r}}}} 

+4224{r,{r,{5,{5,r}}}} 
+896 {r,{r,{r, {5, r}}}} y 




Force 
Gradient 
#2 


6 6 2 f 

48StS-St3 {S%ST}} 
X e 72 

X 


St'* 
155520 


r+5 

^ 41 {5, {5, {5, {5, r}}}} ^ 
+36{{5,r},{5,{5,r}}} 

+72 {{5, r},{r, {5, r}}} 

+84{r,{5,{5,{5,r}}}} 

+126 {r,{r, {5, {5, r}}}} 
V +54 {r,{r,{r, {5, r}}}} / 





Table 3: A collection of integrators with the leading terms in their exactly conserved shadow Hamiltonians. 
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